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Abstract. Let / be a fixed (holomorphic or Maass) modular cusp form, 
with L-function L(f,s). We describe an algorithm that computes the value 
L(f, 1/2 + iT) to any specified precision in time 0(1 + \T\ 7 / 6 ). 



Contents 

1. Introduction 1 

2. Preliminaries 4 

3. A 'geometric' approximate functional equation 7 

4. Proof of theorem [5] 11 

5. Proof of theorem Q] 13 

6. Lemmas 14^31 1^3] and computing c x , y 

7. Numerical integration for analytic functions 2C 
Appendix A. Special functions and lemma I3~2l 21 
References 23 



1. Introduction 

In this paper, we consider the problem of finding a fast algorithm to compute 
the .L-function of a (holomorphic or Maass) cusp form at a point on the critical 
line. We also discuss the problem of computing the large Fourier coefficients of / 
at a cusp. We use very similar algorithms to solve both problems. 

Let r be a lattice in SL(2,R). Let / be a cusp form on r\H. Specifying T and / 
requires (a priori) an infinite amount of data. Therefore we discuss what we mean 
computationally when we say 'given / and a cofinitc volume subgroup T of SL(2, K)' 
in 1231 

The main results can be summarized into the following two theorems: 

Theorem 1. Let f be a holomorphic or Maass cusp form for a congruence subgroup 
of SL(2,Z), and let 7, e be given fixed positive reals. Then for any T > ; we 
can compute L(f, 1/2 + iT) up to an error of 0(T -7 ) in time 0{\ + T» +e ) using 
numbers of O(logT) bits. The constants involved in O are polynomials in and 
are independent of T . 

Theorem 2. Given a lattice T ^ SL(2, R) containing ^ J | j, a (holomorphic 

or Maass) cusp form f on r\H, positive reals 7, e, and a positive integer T, the T th 
Fourier coefficient of f at any cusp can be computed up to an error of 0(T^ 7 ) in 
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time 0(T 7 / 8+e ) using numbers of O(logT) bits. The constants involved in O are 
polynomials in -^jp and are independent of T . 

Our algorithm in theorem [1] is considerably faster than the current approximate 
functional equation based algorithms which need the time complexity 0(T 1+0 ^). 
Our method uses a form of "geometric approximate functional equation" (see (|3. 10|) 
and p.lip ). It can also be used to compute derivatives of L(/, 1/2 + iT) because 
similar "geometric approximate functional equations" can be obtained in these 
cases. 

1.1. Historical background and applications. The problem of "computing" 
values of the zeta function effectively goes as far back as Riemann. Riemann 
used the Riemann Sicgcl formula to compute values of the zeta function and verify 
the Riemann hypothesis for first few zeroes. The Riemann Siegel formula writes 
£(1/2 + iT) as a main sum of length 0{T X / 2 ) plus a small easily "computable" 
error. The first practical improvement for this method was given by Odlyzko 
and Schonhage in [13]. This algorithm can be used to evaluate OiT 1 / 2 ) values 
of C(l/2 + ii) for t G [T, T + T 1 / 2 ], to an error of 0(T~~<) each , using 0(7^/2+0(1)) 
arithmetic operations on numbers of O(logT) bits and using 0(T 1 / 2+0 ( 1 ^ 1 ) storage. 
This method is extremely useful for computing large values of ^"(1/2 + iT) but it 
does not improve the time required to evaluate a single value of £ function. 

Schonhage gave an algorithm to compute £(1/2 + iT) in 0(T 3 / 8+o(1 )) in [17] . 
Heath-Brown improved the exponent to 1/3. Hiary gave an explicit algorithm for 
the exponent 1 /3 using rapid computation of the truncated theta sums in [9] . The 
exponent was later improved by Hiary himself to 4/13 in [S]. The main idea behind 
these algorithms is to start with the main sum in the Riemann Sicgcl formula, cut 
the sum into "large" number of subsums having "large enough" lengths and try 
to compute these subsums rapidly. Fast numerical evaluation of ^(1/2 + it) was 
also considered by Turing [21], Berry and Keating [2], Rubinstein [14] and Arias 
De Reyna [1] et al. 

In the case of higher rank L-functions, the analogue of the Riemann Siegel for- 
mula is given by the approximate functional equation. A detailed description of 
"the approximate functional equation" based algorithms is given by Rubinstein in 
[14] . Booker gave an Odlyzko and Schonhage type algorithm to compute multiple 
valuations of £(1/2 + iT) in [4]. However, algorithms based on the approximate 
functional equation remain the fastest way to compute a single value of L(l/2 + iT), 
when T is large. In the case of L-function associated to a modular (holomorphic or 
Maass) form, these algorithms have 0(T 1+0( ^ 1 ') time complexity. 

Our algorithm in theorem Q] is the first known improvement of the approximate 
functional equation based algorithms in the GL2 setting. 

Computing values of L-functions on the critical line has various applications 
in number theory. It can be used to verify the Generalized Riemann Hypothesis 
numerically. It has also been used to connect the distribution of values of L- 
functions on the critical line to the distributions of eigenvalues of unitary random 
matrices via the recent random matrix theory conjectures. 

The problem of computing L-functions is closely related to the problem of finding 
subconvexity bounds for the L-functions. More generally, improving the "square 
root of analytic conductor bounds" coming from the approximate functional equa- 
tion is of great interest to analytic number theorists. 
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Our algorithm starts with writing L(f, 1/2 + iT) as an integral over the cycle 
{(—1 + i/T)t}. The equidistribution of this cycle in SL(2,Z)\H is closely related 
with the moment J_ T |C(l/2 + it)\ 4 dt (see [16]) and with subconvexity bounds for 
the GL2 L-functions. The results in this paper are closely related to the work of 
Venkatesh [22], where he writes special values of L-functions as period integrals 
over certain cycles. He then uses the equidistribution properties of these cycles to 
get subconvexity bounds for the L-functions. We however use the slow divergence 
of the cycle {(— \ + i/T)t] to get a fast way of numerically computing L(/, 1/2 + iT), 
up to any given precision. 

In the present paper, we have only considered the GL2 case. It will be of great 
interest to generalize the method in this paper to higher rank L-functions . A 
similar (to I3.11j) form of geometric approximate functional equation has been used 
by Sarnak in [16| section 3], in the number field setting. Our method thus could 
also generalize for L-functions corresponding to the Maass forms for SL(2,o), where 
is the ring of integers for a number field K . A more interesting problem will be 
to generalize our technique in the GL„ setting, where the subconvexity bounds are 
also not known. 

1.2. Outline of the proof. We use a geometric method: a suitable integral rep- 
resentation of the L-function reduces the problem to computing an integral of / (or 
rather, /, the lift of / to T\SL(2, M) defined bv !2.2l ) over an approximate horocycle 
of length T. 

Wc then cut the integral over the approximate horocycle into a sum of integrals 
over smaller segments of equal length M . We sort the starting points of the segments 
into sets Si, S-ji ••• such that the points in each set are very close to each other. The 
key point is that for each i, we can compute all the integrals over all segments with 
starting points in Si "in parallel." 

To accomplish this, we use the critical fact that if two points Xi, x% in T\SL(2, R) 
are sufficiently close to one other, then the images Xi(t),X2(t) of these points under 
the time t horocycle flow stay very close to each other "for a long time." The 
approximate horocycles also have a similar property. We quantify these in lemmas 
I4.2l and r5.2l latcr. This is a very special property of the horocycle flow. The geodesic 
flow, for example, does not satisfy such a property 

We return to the description of how to compute the integrals over various seg- 
ments "in parallel." Consider Si, for instance. Let / be a fixed segment with 
starting point in Si and let J be any other segment with starting point in Si . Let 
F(J) denote the value of the integral over the segment J. We use lemmas [4. 2\ 15.21 
to compute / at each point on J using a power series expansion around a corre- 
sponding point on I. In the case of L-value computations, this converts F(J) into 
sum: 

(1.1) F(J)=J2^(J) J rf,,., Hull 

Here, fij,i(t) are smooth functions independent of J and aij are precomputable 
constants depending on J and the sum over i and j is over a set of size 0(1). Wc 
can compute each integral on right hand side of (|1.1|) using 0(M 1+ °( 1 )) operations. 
Therefore, by prccomputing certain constants, we are able to compute all the J- 
integrals in at most a further <3(|Si|) time. 
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This grouping leads to a speed-up by a factor that is roughly the size of the 
groups \Si\. In practice, we cannot make these groups very large, because the 
amount of time for which the points x\{t) and X2(t) stay "very close" to each other 
depends on how close the points x\ and xi were. Hence, if we try to make the 
groups "too large" , the "admissible" length of the segments decreases, resulting in 
an increase in the number of segments. This increases the running time. 

1.3. Model of computation. In practice, specifying a real number completely 
requires an infinite amount of data. Hence to keep the algorithms simpler and 
clearer, we will use the real number (infinite precision) model of computation that 
uses real numbers with error free arithmetic having cost as unit cost per operation. 
An operation here means addition, subtraction, division, multiplication, evaluation 
of logarithm (of a complex number z such that | arg(z) | < 7r) and exponential of a 
complex number. 

Our algorithm will work if we work with numbers specified by O(logT) bits. 
This will at most add a power of log T in the time complexity of the algorithm. We 
refer the readers to [23J Chapter 7] for a brief discussion of the floating point error 
analysis for this algorithm if we use numbers specified by O(logT). We refer the 
readers to [19} Chapter 8] and [20] for more details about the real number model 
of computation. 

1.4. Outline of the paper. A brief account of the notations used in this paper is 
given in section [5] The algorithm uses a type of "geometric approximate functional 
equation" . It is discussed in detail in section [3l The lemma 13.21 is the main tool 
used in deriving the "geometric approximate functional equation" for the Maass 
form case. We will give a proof of it in appendix [A] We will give the proof of 
theorems [2] and [1] in sections [4] and [5] respectively. We will discuss the proof of the 
"special property" of the horocycle (quantified in lemmas 14.21 15.21) and using it to 
simplify the integral (quantified in lemmas |4.3H5.3|) in section[6]in detail. In section 
[7J we will give a simple algorithm for numerical integration of "well behaved" real 
analytic functions. 
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me this problem, as well as for all the help and support he gave me. I learnt a 
great deal about this subject from him. I would also like to thank D. Bump, K. 
Soundararajan, M. Rubinstein, P. Michel, A. Booker and A. Strombergsson for the 
encouragement as well as for helpful discussions. I like to thank F. Brumley, G. 
Hiary and H. Then for providing helpful references. I also like to thank M. Raum 
and A. Parthasarathy for helpful editorial comments. 

Majority of the work in the paper was carried out during my PhD at the Courant 
Institute and Stanford University, and also during my stay EPFL and MPIM, Bonn. 
I am very thankful to these institutions for allowing me to visit and carry out my 
research. 



2.1. General notation. Throughout, let T be a lattice in SL(2, Z), unless specified 
otherwise. In the proof theorem [21 T will denote a lattice in SL(2,K) containing 
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Throughout, let T be a positive real. In section |4j however T will be a positive 
integer. We will denote the set of non negative integers by Z+ and the set of 
non-negative real numbers by R+. 

We will use the symbol <C as is standard in analytic number theory: namely, 
A <C B means that there exists a positive constant c such that A < cB. These 
constants will always be independent of the choice of T . 

We will use the following special matrices in SL(2,R) throughout the paper: 

, 01 , ... ( 1 t \ , . / evl 2 \ / cos<9 sin6» \ 

e(x) will be used to denote exp(27rix). 

Let / be a cusp (Maass or holomorphic) form of weight k on r\H. The weight 
corresponding to a Maass form will be 0. Let us define a lift / of / to r\SL(2,R) 
by / : T\SL(2, R) — > C such that 

(2-2) f(( : b J ))=(ci + d)- k f 



c d 

Let x be any element of SL(2, R). We will frequently abuse the notation to treat 
x as an element of r\SL(2,R). i.e. we will often denote the coset Tx simply by x. 

2.2. Real analytic functions on T\SL(2, R). Let x be an element of SL(2, R) and 
let g be a function on T\SL(2, R), a priori g(x) does not make sense but throughout 
we abuse the notation to define 

g(x) = g{Tx). 

i.e. g(x) simply denotes the value of g at the coset corresponding to x. 
Let </> be the Iwasawa decomposition given by 

<f> : (t, y, 6) e R X R X R -> n(t)a(y)K(9). 

Recall that <j> restricted to the set R x R x (— ir,ir] gives a bijection with SL(2,R). 

Definition 2.1. Given e > 0, let il £ = (— e, e) x (— e, e) x (— e, e) and U e = </>(il e ) C 
SL(2,R). 

Let us define the following notion of "derivatives" for smooth functions on 
r\SL(2,R): 

Definition 2.2. Let g be a function on SL(2,R) and x any point in SL(2,R). We 
define (wherever R.H.S. makes sense) 

d d 
q^9(x) = g^\t=09(xn(t)); 

d d 
^-SW = -Q t |t=o g(xa(t)); 

£-g{x) - ^ | t=0 g(xK(t)). 

Sometimes, we will also use di to denote 
Given /3 = 02, ffs), let us define d l3 g(x) by 

flPl f)P2 

(2.3) d p g{x) = B a-g{x). 

K ' yK ' dxf 1 dx 2 2 dx§ 3 v ; 
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For /3 as above we will define 
and 

1/31 = \Pl\ + \fa\ + \P3\- 

We now define the notion of real analyticity as follows: 

A function g on r\SL(2,R) will be called real analytic, if given any point x in 
r\SL(2,R), there exists a positive real number r x such that g has a power series 
expansion given by 



(2.4) g(xn(t)a(y)K(9)) = £ ^M t ^ y ^ 

for every (t,y,9) G U rx . 

Let us use the following notation for the power series expansion. 



Definition 2.3. Let y, x G SL(2, R) and t, y, 9 be such that y = xn(t)a(y)K(9) and 
{J3 1 ,/3 2 ,l3 a )=j3€l? + define 

(y-xf = t h y^Q h . 
Hence we can rewrite the Equation (|2.4j) as 

9 (y)= E ^(y-xf- 

Throughout, we will assume that for a cusp form /, all the derivatives of the lift 
/ (of /) are bounded uniformly on T\SL(2,R) by 1. In general it can be proved 
that given a cusp form /, there exists R such that ||9^/||oo "C R^', see [221 section 
8.2]. The case when R > 1 can be dealt with analogously. The assumption that all 
derivatives are bounded by 1, allows the proofs to be marginally simpler. 



2.3. Input for the algorithm. As mentioned before, to specify T and / com- 
pletely, a priori we need an infinite amount of data. Throughout the paper we will 
assume that T is given in terms of a finite set of generators. Similarly we assume 
that modular(holomorphic or Maass) cusp form / on r\H is given by first O(logT) 
of its fourier coefficients. 



2.4. Computing values of /. In practice, values of / and its derivatives can be 
computed very rapidly. See [3) for a method for effectively computing /. It can be 
easily shown that given any x and any fixed 7, one can compute d" f(x) up to the 
error 0(T~ 7 ) in 0(T°^) time. Here the constant involved in O is a polynomial in 
|/3 1 and 7. In this algorithm we only compute values of d@f(x) for |/3| <C 1. Allowing 
0(T , °( 1 )) time for each valuation of / does not change the time complexity of the 
algorithm. See [23l chapter 7] for explicit details about it. Hence for simplicity we 
will assume that each value of / (or /) or any of it's derivative can be computed 
exactly in time 0(1). 
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3. A 'GEOMETRIC' APPROXIMATE FUNCTIONAL EQUATION 

Most existing algorithms to compute a general L-function start with an approx- 
imate functional equation. This is a generalization of the Riemann- Siegel formula 
for a general L-function. We refer the readers to P31 Section 3] for a detailed 
discussion of the approximate functional equation based algorithms. 

Our algorithm however starts with a 'geometric approximate functional equa- 
tion'. Which means that we write L(f, 1/2 + iT) as an integral over a 'nice' curve 
of hyperbolic length 0(T 1+0 ^) up to any given polynomial error. This corresponds 
to the usual discrete sum form of the approximate functional equation. The deriva- 
tion of the approximate functional equations used here have geometric motivation 
behind them. Hence we call these equations as "geometric approximate functional 
equations". The 'geometric approximate functional equation' for _L(/, s) is given 
(for holomorphic and Maass forms respectively) by p,10|) and p. lip . 

Let r be a lattice of SL(2, Z) and let / be a modular (holomorphic or Maass) 
cusp form of weight k. In the Maass form case, let / be an eigcnfunction of the 
hyperbolic laplacian A = —y 2 (d 2 + d 2 ) with eigenvalue 1/4 + r 2 on r\H. Here r is 
any positive real number or it is a purely imaginary number with \r\ < 1/2. 

We will further assume that / is an eigenfunction of some Fricke involution, as 
is the case for all the newforms. This will quantify the exponential decay at zero. 
This implies that there exists a positive constant C\ and a real constant C2 such 
that: 

(3.1) f { -9l)=c 2 z k f{z). 

In this paper, for simplicity let us assume that / is either holomorphic or an even 
Maass form. The algorithms will be analogous for the odd Maass cusp form case. 
For an even Maass cusp form, we will use the following power series expansion : 

(3.2) f(z) = ^f(n)W r (nz). 

n>0 

Here W r (x + iy) = 2- s /yKi r (2ny) cos(27ra). Here K ir denotes the /T-Bessel function. 
See [10] for more details. The explicit Fourier expansion for the holomorphic cusp 
forms is given by 

(3.3) f(z) = J2f(n)e(nz). 

n>0 

Given an (even Maass or holomorphic) cusp form / of weight k, let us define the L 
function by 

/(") 



Definition 3.1. L(f,s) = £~ 1 nS+ T> 1)/2 . 

The usual integral representation for L(f,s), corresponding to a holomorphic 
cusp form is given by 

/>oo 

(3.4) / f{it)t s+i - k -W 2 dt = (2 7 r)-' , -( fc - 1 )/ 2 L(/, s)T (s + (k - l)/2) . 



Using the condition (|3.1j) and the exponential decay of / at 00, the integral on the 
left hand side of (|3.4j) converges for all s in C. Hence we can use (|3.4j) to compute 
the values of the L functions on the critical line. However, the gamma factor on 
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the right hand side of (|3.4[) is decaying exponentially with T (here T = Im(s)), 
therefore we need to compute the integral to a very high precision in order to get 
moderate accuracy in the computations of the L-valucs. This problem has been 
tackled before in several ways. A solution to this has been suggested in [IT| and 
worked out by Rubinstein in [15] Chapter 3] . 

We will proceed in a different fashion. We change the contour of integration 
that will add an exponential factor which will take care of the exponential decay of 
the gamma function. We will derive a 'geometric approximate functional equation' 
using this idea in this section. 



3.1. A geometric approximate functional equation for L(f, s). Let s = 1/2+ 
iT and a = — 1 + 

As stated before, we change the contour of integration in (|3.4p so as to add 
an exponential factor which will take care of the exponential decay of the gamma 
function. The case of holomorphic modular forms is easier to deal with, hence we 
will consider it first. 

Holomorphic case. : Let / be a holomorphic cusp form of weight k > 0. 



(3.5) 

f{at)t s+ ^/ 2 - l dt 



V/(n) / exp(2TTinat)t s+{k - 1 ^ 2 - 1 dt; 

-, JO 



nS +(fe-l)/2 (_ m ) S +(fc-l)/2 J z= _ iaR+ 

~ [Z7T) 2- T)S +(fc-i)/2 (_ ai ) s +(fc-i)/2 / „ CX P^ z > z az 



n=0 



= (2.)-^-^ ( Jffi£ 1)/2 r(« + (* - D/2). 

Here, the branch cut for the logarithm is taken along the negative real axis. 
Hence the argument of the logarithm takes values between —tt and n. With this 
choice of the logarithm, the factor ^_ ai yl( k -D/2 grows like 0{e" T / 2 ) as T — > oo, to 
compensate for the exponential decay of T{s + (k — l)/2). 

Case of Maass forms. : We use a similar idea for the non-holomorphic case, let 



We will choose a suitable value for T\ later. 



RAPID COMPUTATION OF L-FUNCTIONS FOR MODULAR FORMS 



9 



/ f( ai t)t'- 3 ' 2 dy = /(«) / Wrinaxty-^dt; 

J ° n>0 J ° 

=E^r^((- i +^-)^- 3/2 ^ 

Tl>0 J0 1 

1 - S " 1/2 -A, 
- 2T ' E ^T72 ^ cos(Tit)K lr (t)t s ~ 1 dt 



cos(2Trt)K lr (2Tr—)t s - 1 dt; 



We summarize the result as: 



(2tt) s 

v ; n>0 



c y r j^^ 1/2 



(3.6) / f( ai t)t s -^ 2 dy = 1 L{S, a) / cos(T 1 t)A^ r (t)f ;T " 1 / 2 ^. 



n 



In lemma [3721 we will show that we can choose T\ = 0(T) such that the integral 
on the right hand side of p.6[) is not too small. 



Lemma 3.2. Let r be any fixed complex number such that \Im(r)\ < 1/2 and 
T,Ti > 0. T/ien, 

(3.7) cosiT^K^tyr-^dt = 2 ^-V2 T{ ^±E±l n Z^±E±l ) 

Jo 2 2 

ir + iT + h —ir + IT + h 1 9x 
F( 2-^,3,-2?). 

iiere F denotes the hypergeometric function. Moreover there exists a computable 
constant C = 0(1) depending on r such that for T\ = C'T , 

T( ir+jT+\ , T , -ir + iT+\ ir + iT + \ -jr+jT+j 1 
I 2 ' 1 2 jl 2 ' 2 ' 2' 1 J 



T7ie implied constant is non zero and depends only on r. 



A proof of the lemma 13.21 can be found in [5]. However, to keep the paper 
somewhat self contained we will prove lemma [3~2l in appendix |A"1 It will allow us to 
take T\ = C'T for some C = 0(1). Hence for this choice of Ti, once we compute 
Jo f{otii)t s ~ l dt up to an error of 0(T -7 ), we can compute L(f, s) up to an error 
of OCT 1 "?). 

We next use the exponential decay of / at and at oo to get: 

Lemma 3.3. Given any cusp form f , a real number Sq and positive reals T, 7, if 
a = —l + i/T and s = so + iT then there exists a computable constant c = 0(1 + 7), 
independent of T such that for any <ii, g?2 > c: 

/>oo pdzTXogT 

(3.8) / f{at)t s - 1 dt= / f(at)t s - 1 dt + 0(T-' 1 ). 
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Proof. Using the exponential decay of / at oo, we get that for t > 1, 

\f(at)\ = \f(-t + it/T)\ « exp(-af/T) 
for some positive constant a. Hence for all t > T, we have 

\f(at)t s - l \ < t 80 ' 1 exp(-at/T). 

Hence, 

t* - 1 < exp(at/2T) (s - l)logf < at/2T ^ 2T(s - 1) log i/a < t. 

It is easy to see that there exists a constant c' depending only on a and so such 
that the above condition holds for t > to = c'TlogT. This implies that for t\ > to, 



(3.9) |/ f{at)t s - 1 dt\<^ I cxp(-at/2T)dt = (2T/a) exp(-a*i/2T). 
Jti Jti 

Hence for any ci > 2(c' + 1)(1 + l/a)(l + 7), we get that 

\f(at)t s - l \dt^.O{T-<). 



' Cl TlogT 

Similarly using (|3 . 1 [) , we get that for t < 1 < T, we have 

\f(at)\ = \f(-t + |)| = M^ t -*|/(-Ci/H + |))l « ex P (-Cxa/2Tt). 

We deal with this case exactly as in the previous case to get a suitable constant 
c. □ 

Applying Lemma 13.31 to (|3.5[) we get that given any T, 7 > there exists a 
computable constant c such that for any to < cT ^ T , and for any ii > cTlogT, 

(2 7T )s+(k-l)/2(_ -\s+(k-l)/2 ftt 

(3.10) L(f,s)= W r(5 + ( fe_ 1 } /2) / /M)^- 3 )/ 2 dt + 0(T--). 

and similarly after applying Lemma [3.3l and lemma EOl to (J3TBJ) we get that given any 
T, 7 > there exist computable constants c and C = O(l) such that if T\ = CT 
then for any to < eT ^ T , and for any fi > cT log T, 

(3.11) L(/, .) = ^i^^^ £ /(a!*)*-"/ 31 * + O(T-). 

Here, C(Ti) = J Q °° cos(T 1 t)K ir (t)t s - 1 dt and T x is chosen such that C(T X ) > 
T -1 . Notice that the integrals on the right hand side of (|3.10[) and (|3.11[) can 
be taken over curves of hyperbolic length w TlogT. We call (|3.10[) and (|3 . 11 1) as 
geometric approximate functional equations for L-functions corresponding to the 
holomorphic and the Maass forms respectively. It is easy to see that the integrals 
on the right hand side of the equations (|3.10[) and (|3.1ip can be computed directly 
in 0(T 1+e ) time, given any e > 0. Our algorithms for computing L(f, s) start with 
these "geometric approximate functional equations" . 
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4. Proof of theorem [2] 

Our main aim is to prove theorems Q] and [2] Both the algorithms are very similar 
and they use basically the same idea but the algorithm for theorem [2] is marginally 
simpler. Hence we discuss it first and in the next section we will deal with theorem 

HI 

Let r be a lattice in SL(2, K). Given a real t, let 

ei 
e~i 



«*(*) = 
and 

»(*) = 



1 t 
1 

Let / be a holomorphic cusp form on r\SL(2,R) of weight k. In this section we 
will only discuss the holomorphic case. The algorithm in the case of a Maass form 
/ will be analogous. Recall the lift / : T\SL(2,R) -> C of / defined in ((221) as: 

5 fa b \ . , ,._ fe , / ai + b 



c d J \ci + 

Recall that we have assumed that / is "well behaved", i.e. that / has bounded 
derivatives. The Fourier expansion is given by ()3.3[) : 

n>0 

Hence we get that for any positive integer T, 

e- 27r f(T) = [ f(x + i/T)e{-Tx)dx 
Jo 

(4.1) = f T fe / 2 - 1 /(a(-logT)n(t)) e (-t)di. 

Jo 

Notice that the integral on the right hand side of (|4.1[) is an integral of a well be- 
haved smooth function on a horocycle of length T. Hence a priori we can 'compute' 
it in 0(T 1+ °^) time up to the error at most 0(T~ 7 ), for any given 7. This will 
denote the corresponding (to (|3.10|) and p.lip 1 ) "geometric functional equation" in 
this case. We will now explain a method to compute the T th Fourier coefficient 
faster than O(T). 

Let 77 be any positive number. We will write the integral (|4.ip as a sum of 
integrals over horocycles of length T v each. For simplicity let's assume T v is an 
integer. Hence we have, 



f T IT 1 "".!-! pT v 

(4.2) / f(x n(t))e(-t)dt= V / f{x D n(jT v + t))e(-t)dt 

Jo J=0 Jo 

+ [ f(x n{t))e(-t)dt. 

Here xo = a(— logT). The second integral on the right hand side of (|4.2[) is an 
integral of a smooth well behaved function on a horocycle of length at most T v . 
Hence given 7, e > 0, we can compute it up to an error of 0(T f ) in time 0(T v+e ), 
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using proposition [73] In practice, e will be a fixed "small" real number and 77 will 
eventually be chosen to be 1/8. 
Let 

M = T' n . 

Let us define Ii(x, /) by 

Definition 4.1. Given smooth function g on r\SL(2,R), x in r\SL(2,R) and a 
non negative integer I we define 

(4.3) Ii(x,g)= t l g(xn(t))e(-t)dt. 

Jo 

Hence we can rewrite the equation (|4.2[) as 

/ f(x Q n(t))e(-t)dt = V h{ XQ n{jT>),f) 

+ / f{x n{t))e{-t)dt. 

As mentioned in subsection ll.2l the slow divergence of the horocycle flow will imply 
that a lot of the pieces (in (|4.2|1 ) of the horocyle will be very close to each other. 
We will group the pieces, which stay "very close" to each other and try to compute 
the integrals on all the pieces in each group "in parallel" . 

The slow divergence property of the horocycle flow used here, is quantified in 
the following lemma 1431 

Lemma 4.2. Given any e > and r\ > 0, for any rj' > 2r) + e and for any x,y 
such that x~ 1 y € U T - V > , we have a constant c, independent of r/ and rj 1 such that 
yn{t) e xn(t)U cT -. for all < t < T' 1 . 

Let e be any positive number. We take points {xon(jM),0 < j < T 1_v } and 
"reduce" the points to an approximate fundamental domain. Then we "sort" the 
reduced points into the sets Si, S2, Sn such that x, y e Sj =>■ x~ x y e U T -(2^+ e ) . 
For each j, let us choose a representative Vj from each Sj. 

It is easy to sec that N < T 6r ' +3e logT. 

If y € Si, then lemma l4~2l implies that the points yn(t) and Vin(t) stay "very 
close" to each other for all < t < T v . Therefore we use the power series expansion 
around Vin(t) to compute values of / at yn(t). Hence, we get the following lemma: 

Lemma 4.3. Given 7 > 0, e > 0, any 77 > and x, y <E Si for some i, then we have 
constants c XtV jjj and d such that 

(4.4) io(y,/)= E E^r J,(a? '^ /)+0(T ~' Y) - 

\p\<d,p&\ i=a 

Here, d = 0((1 +7)/e). The constant involved in (|4.4|) is polynomial in d. 

Notice that the equation (|4.4[) is the explicit form of (jl.ip in this case. Let us 
observe that the integrals involved on the right hand side of (|4.4[) depend only on x. 
This lemma will allow us to compute the sum Io(y, f) in parallel for all the points 
y in Si in 0(\S. t \ + M 1+£ ) steps. 
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Lemma l4~3l implies that we can get a number d = 0((1 + 7)/e) such that 

.T N d 

/ f(x on (t))e(-t)dt = J2 E T. c -^ir Ll ^ dfi ~f) 

m=ly£S m \p\<d,pez% 1=0 

(4.5) + f j(x o n{t))e{-t)dt + 0{T-~<). 

We will prove the lemmas H . 2 1 and |4~31 in section [5] Let us complete the proof of 
theorem [5] using lemma 14.31 

Proof of theorem [H Notice that (|4.5[) , along with (|4. 1[) give us a method of com- 
puting f(T) up to 0(T -7 ) error. Let us compute the time spent in this algorithm. 

It is easy to see that "reducing" each Xi to an approximate fundamental domain, 
can be done in 0(log T) time. There are many standard reduction algorithms 
available. We refer the readers to [23] for a form of the reduction algorithm. The 
whole "reducing" and sorting process has also been discussed in detail in |23[ chapter 
7]. 

As mentioned earlier, N, the number of sets {Si} is « 0(T 6ri+3<i logT). Hence 

the total time needed in reducing all the x;'s to the fundamental domain, then 

sorting them into sets Sj's and picking a representative Vj from each Sj requires 
^ T i-v + T ev+3c^ logT ) stepg ^ 

Notice that in equation (|4.5[) . the integrals Ii(v m ,d" f) arc independent of the 
choice of y. For each y € S m , there are 0(1) terms in the right hand side of equation 
(|4.4[) . In section 16.11 we will show that for fixed y and v m , we can compute each 
c v m ,y,/3,l's in O(l) time. Hence we can precompute the constants Cy y pj, for all y 
in S mi and for all m in 0{T l ~ v ) steps. The constants involved are polynomial in 

(1+7)A- 

Recall that M = T n . For a fixed v m , computing (v m , f) for all |/3| , I < d takes 
0(T v+e ) operations (using proposition 17. ip . Hence, using Lemma (|4.3[) . we need 
|S m | more operations to compute Io(y, /) for all y € S m . The maximum number of 
the sets S' m s is 0(T 6,,+3e logT). Hence the total time required to compute Io(xi, f ) 
for all Xi is (up to a polynomial factor in (1 + j)/e) at most 

N 

(T n+e + \S m \) < T" +£ T 6 " +3e + T 1 -'". 

m— 1 

Notice that the extra log factors can be absorbed at the end in the exponent e. 

The second integral on the right hand side of (|4.2|) . can be computed up to an 
error of at most (9(T -7 ), using at most 0(T T ' +e ) operations. Hence the total number 
of operations needed to compute J f(x n(t))e{t)dt up to an error of 0(T -7 ) is 

(4.6) (( r ^ T 6r,+3<r + 

The optimal value for rj is 1/8. This implies that given any e, 7 > and real so, 
f(T) can be computed up to an error at most 0(T -7 ) using at most 0(T 7 / 8+4c ) 
operations. The constant involved is a polynomial in (1 + j)/e. □ 

5. Proof of theorem [1] 

The proof of this theorem uses an idea very similar to the proof of theorem 
[U Let T be a lattice of SL(2,Z). Let / be a (holomorphic or Maass) cusp form 
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on r\H. Let a = -1 + i/T and c > be the constant in ([3~10| (or in (f3~TT|) for 
the Maass form case). Let s = sq + iT and to = cTlogT. Let so be some fixed 
real number. Using the "geometric approximate functional equations" (|3.10[) and 
(|3.11[) , it is enough to give an algorithm to compute 

r*i 

f(a£)t 8 - l dt, 



for some suitable ti > cTlogT. Recall that in the Maass form case, we need 
an algorithm to compute f(ait)t s ~ 1 dt. However, the algorithm for comput- 
ing J t 1 f(at)t s ~ 1 dt can be trivially generalized to an algorithm for computing 
J t * J f(ait)t s ~ 1 dt. Therefore, in this section we will only give an algorithm to com- 
pute f£ /(at)< s_1 dt. 

Let T] be any non-negative number. We proceed in a similar manner as in the 
previous section and write the above integral as a sum of integrals over segments of 
hyperbolic length sa T v . In order to do so, let us first define the following quantities: 

Definition 5.1. 

2/o = ha, 

bj = (i + r- 1+ yt , 
yj = M- 

Let M x be an integer such that M x = 0(T 1 -'» logT) and 

&Afi < cT log T < b Ml +i 

Let 

tl = &Mi+l- 

Here yo denotes the starting point of the approximate horocycle. The points 
Ui denote the starting points of segments. The hyperbolic distance of the segment 
between z/, and j/j+i is w T v . Mi is the total number of segments required. 

Recall that c is chosen such that it satisfies conditions in (|3.10j) / p. lip . Hence it 
is enough to give an algorithm to compute 

i-ti 

f(at)t s - l dt 



up to an error at most 0(T 7 ) 

rti j^i r b s T- 1+n 

f(a(bj+t))(bj+ty- l dt 



/•tl M l r- 

/ fiaty-Ht = £ / 

J to J=0 Jo 



Mi „T _1 + * 7 

(5.1) =E^/ /(a(6,(l+t))(l + t) s - 1 & 
Let us change the variable to u = tT. Hence we can rewrite (|5.1|) as 

(5.2) / f(at)t s - 1 dt = T- 1 Y J b° / 1(^^(1 + u/T))(l + u/T) s - 1 du. 
Jt j=0 Jo 

Here M = T v . Let 

(t/T)i -(tT)i 



(5.3) n(t) 



(t/TY 
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n(t) denotes a lift of the curve {at} to SL(2,R). Notice that 

fu -w uu m x ( (l + u/T) 1 / 2 -u(l + u/T)" l /2 
+ = I o (l + u/T) 1/2 

and that it is independent of j . Let 

(5-4) «(«)=^ Q (i + u/T )-V2 J- 

Let = be a lift of the point yj to SL(2,R). We rewrite equation (|5.2[) as 

r-ti 



/(at)^" 1 ^ 

to 

Mi <>A/ 

T -1 bj / /(a(6 J (l + u/T))(l + w/r) s " 1 rfu 

Mi „T" 

T fc/2-i 6 »-*/a / (1 + u /T) s - fc / 2 - 1 /(/ t (6, (l + u/T)))<ft 

Ml nT* 7 

T */2-i^ 6 '-*/a / (i + u /T)- fc / 2 - 1 /(«(6 i )«(6 J -)- 1 «(6i(l + «/T)))dt 

.•_n JO 



J'=0 

(5.5) 

Mi 

= T fc /2-i^ 6 -^ jLo(/ >,) 
i=o 

Here given a smooth function 5, we define 

Li(g,x)= / u l g{xu(u))(l + u/T) s - k / 2 - 1 du 



We make use of the fact that we are integrating on an "approximate horocycle" . 
In particular, that it has slow divergence. We quantify this result in the following 
proposition. This lemma is analogous to lemma [ 



Lemma 5.2. Let 77 > 0,e > be such that e < 1 — 3r/ ; and x,y G SL(2,R) such 
that x~ x y G Ut-^-c , then we have a constant c' such that 

yu{u) G xui(u)U c it-c for all < u < T v . 

c' can be chosen independent of rj and T. 

The lemma [5T21 will be proved in section [BJ Let's 'reduce and sort' the points 
{xj} into N = 0(T 6t '+ 3£ logT) groups Si, S N such that 

x,y G Si =>■ x~ 1 y G Ut~2v->- 

It is again easy to see that the number of the sets Si's is 0(T 6r?+3£ logT). Let 
us also choose fixed representatives from each Si. For each group Si, let us 
try to compute the inner integrals (in (|5.5[) ) corresponding to the points in Si "in 
parallel" . 

Let x, y G Si, then we use power series expansion around points xu(t) to compute 
f(yoj(t)). Hence we get, 



16 



PANKAJ VISHE 



Lemma 5.3. Given any e,t], 7 > 0, such that e < 1 — 3r), and 1,1/6 Si for some 
integer i, then there exist constants e x , y ,p,i and d! , independent of T , such that we 
can write 

Lo(f,y)= E E^^(^/» + o(t-^). 

\P\<d>,l3<EZ 3 + 1=0 

Here d' <C (1 + 7)/e. The constants involved in O are polynomial in ((1 + j)/e). 

Lemmas 15.21 and 15.31 will be proved in section [6l Lemma 15.31 gives the explicit 
form of (jl.lj) in this context. 

Lemma 15.31 implies that there exists a computable constant d! — 0((1 + 7)/e) 
such that 

(5.6) 

tl N d' 

/ f(at)t°-'dt=J2 E a v E E £! ^ L ^/> m ) + 0(T^) 

Here, for each y, let y = x n for some n then the constants a y are defined by 
a y = T k / 2 ~ 1 bn k ^ 2 ■ The lemma 15.31 and (|5.6p convert the problem of L- value 
computation into the problem of computing Li(d@ f ,Vi) for each Vi, and for all 
I, \0\ < 1. Recall that 

L l (f,x)= / M '/> W (u))(l + U /r) s - fc/2 - 1 ^. 



Notice that each integral L/ (d 13 f, v{) is an integral of 87 on a segment of hyperbolic 
length sa T'', hence we can compute it up to an error of 0(T -7 ) in time 0(T V ). 
More rigorously, we prove: 

Lemma 5.4. Given an integer I > 0, a vector f3 £ Z^_, x G r\SL(2,M) and given 
positive reals 7 and e, we can compute Li(d^f,x) up to an error at most 0(T / ) 
in 0((T' ) ) 1+£ ) time. Here, the constant in O is a polynomial in I, \j3\ and (1 +j)/e. 

Proof. Let gi be the function defined by gi(u) = u l (l + u/T) s ~ k / 2 ~ 1 . It is easy to 
see that there exists a constant C, independent of I such that for any n G N, and 
for all < u < T v } 

(5.7) \d n gi (u)\^n\{lC) n {l + u l ). 

If we prove that for any n G N, fixed (3 G and for any x G T\SL(2,K), there 
exists a constant D independent of x and I such that, 

(5.8) \d n { gi (u)dPf(xu(u)))\ « n\(lD) n {\ + «') 

then using a proposition 17. 11 we get the result. 
Recall that 

{l + u/Tf/ 2 -uil + u/T)- 1 ' 2 
{l+u/T)- 1 ' 2 

(5.9) =n(-u)o(log(l + «/T)). 
Let < uq < T v . It is easy to verify that 

uj(u +t)= uj(u Q )n{-t/(l + u /T))o(log(l + (t*o + t)/T) - log(l + u /T)). 



ui(u) 
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Recall that we have assumed that for any j3 in and any x € r\SL(2, 
dPfix) < 1. Therefore, we can use the power series expansion for / (see 
to get 

\0'\=o,/3ez 2 + xo 

(5.10) (-t/(l + u /T) f (log(l + («o + t)/T) - log(l + Uo/T)f* . 

Let us differentiate n times with respect to t. Notice that if j3[ + f3' 2 > n, then 
0»|t=o(-t/(l + u /T))^(log(l + (uo + t)/T) - log(l + u /T))& = 0. Hence the 
only contributions come from the terms for which |/3'| < n. We use this fact and 
well behavedness of / to get 



d n \ u = Uo d f(xu;(u/T))\ « K< 

\P'\<n 



n 3 . 



We use this bound along with (|5.7[) to get (|5 ,8[) and hence the lemma. □ 

Proof of theorem [IJ Using (|5.6p , lemmas 15.31 and 15.41 and following exactly similar 
steps as in the proof of theorem^ we get the result. □ 

6. Lemmas 14.31 15.31 and computing c x ^ y ^ t i 



Lemmas 14.21 and 15.21 are analogous. Similarly lemmas [4.31 and 15731 arc analogous. 
Our main goal is to prove them in this section. 

Let x,y be any two arbitrary points in L\SL(2,R) and let y — xA where A 



r s 



. Given any 5, let Us be the S neighborhood ball around identity defined 



We will use the following easily provable fact to prove the lemmas 14.21 and 15.21 

Fact 6.1. There exist absolute constants C\ and C'2 such that for any < S < C\, 
and any A in SL(2,R) such that \\A — /||oo < o", we have that A € Uc 2 8- 

Let us now start with the proof of the lemma 



Proof of lemma \J7^ The proposition is equivalent to proving that [n{— t)An(t)) £ 
U c t-> for some c. Using fact 16.11 it is enough to prove that \\n{—t)An(t) — /||oo <C 
T~ e for all < t < T v where ||^T||oo denotes the usual infinity norm. But we have 
that 



3.1) n(~t)An(t) 



p — tr (p — s)t — t 2 r + q 
r tr + s 



It is easy to see that for T sufficiently large, if | \A— 1\ \ QCl < T 2ri e then \\n(—t)An(t)- 
I\\oo < T~ e for all < t < T r > . Hence using ECU we get the result. 

□ 

If we look at the the equation (|6.1[) and compute the N, A, K coordinates for 
n(— t)An(t), we get the following immediate corollary 
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Corollary 6.1. Using the same notation as in the proof of lemma \J^ given any 
j3 in 1? + and any k £ Z+, we have constants c x , y ,i3,k such that 

\c xM \ « 2 fe fc 3 /3!T- fe ("+4) 

and 

(yn(t) - xn(t)f = ^ c x ,j, l/3 , fe i fe . 

This also implies that given any 7, e > 7 we have a constant d = 0((1 +7)/e) such 
that for allO<t< T" 7 , 

, s (yn(t) - xn(t)Y 1 h , ^ 

(6-2) pi = J\ E c *M tk + 0{T-i). 

Proof. Given any ^ ^ ^ in SL(2,K), we have 

(6.3) ^=n(^±fM-log(r 2 + S 2 ))A-(tan- 1 (-^ )) . 

We use ()6.3|) and equation (|6.1[) to get 



n (-t)A»(f) =n{ r(P-tr) + (tr + S )(( P - S )t-t*r + q ) 

r 2 + (tr + s) 2 



(6.4) a (_iog( r 2 + (t r + s )2) )ir ( tail -i(_^_)). 

Let fc llft>X)1 ,(t) = ( KP-tr) + (tr+j((p-,) f -tV +g) ) ^ 1; h2Axy{t) = ( _ log(r 2 + (tf + 

s) 2 ))' 32 , and h 3i B 3 , x ,y(t) = (tan"^--^))' 33 . Hence we have 

(yn(i) - xn(t)) fi = hi t p ltXt y(t)h2,0 2 , x , y (t)h 3 ,p s , x , y (t). 

As /IU < r-f 2 ^), we have that \r\,\p-s\, \l - p\, |1 - s|, M < T-( 2l ?+ e ). 
Using these bounds, it is clear that 

(6.5) h< CkxJ°) « /3i!"!2"T-(" +e / 2 )". 

Similar results hold for h2^ 2 , x ,y and ^3,^3, Taylor series for hi^ 2tXt y, ft.2,^ 2 ,z,y 
and /i3^ 3l a:,i/ gives us that 

rfifn r v /4i\^(o)^L,,(o<L,,( ) 

(6.6) C XtVtPi i - an ■ 

/3>£Z%,\P'\=1 

This gives us an algorithm to compute c x>y> B,i- Use (|6.5[) to see that 

c x , y , B ,i « 2^ 3 /3!T-(" +£ / 2 )'. 

Using this bound for c x ,y,p,i we can prove that the radius of convergence for the 
power series X]fcez + c ®,v,pM k ' 1S a * l eas t 0{T e / 2+ri ). Hence using the real analyticity 
of h j>0uX>y for j = 1, 2, 3, we get that for all < t < T''+ e / 2 , we have 

(yn(i) - xn(t)) 13 = ^ c x , y ,B,it l . 



Using this equation and the bounds on c xv g 1, we get (|6.2p . □ 
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Proof of lemma \5.SX Recall 
Notice that 

uj{u) = n{—u) + Err. 
Here H-ErrH^ < T~ 1+2r i for all < u < T r > and 

lq~ x {u) = n(u) + Erri. 

where H-Er^Hoo <C T~ 1+2v for all < t < X" 7 . Using result the proposition is 
equivalent to proving that 

||w- 1 (u)Aa;(n)|| 00 <T- e 

for all < u < T n . However, using the above estimates we get, 

\\uj- 1 {u)Auj{u) - n(u)i4n(-ti)||«> < T- 1+3 ". 

Hence for e < 1 — 3rj, using the same technique as in the proof of lemma we 
get the result. □ 

Again, if we compute the NAK coordinates of uj~ 1 (u)Aui(u), we get following 
corollary for any /3 € 

Corollary 6.2. Using the same notation as in the proof of lemma [S7Sl given any 
j3 in Z\, we have constants e. x ,y,$,k such that 

\e xM \ « 2 fc fc 3 /3!T- fe ("+*) 

and 

(yuj(u) - xuj(u)) p = 2J e.x, y ,p,ku k ■ 
fcez + 

This also implies that given any 7, e > 0, we have a constant d! = 0((1 + j)/e) 
such that for all < u < T' n , 

(yu>(u) — xlo{u))P 



-^X)w, fc u* + 0(T-T). 



K M fe=0 

Proof of lemma \J7S\ Let y = be as defined at the beginning of this section, we 
use power series expansion for / around points xn(t) to compute the value of / 
at points yn(t). Using lemma and the "well-behavedness" of /, we get that 
for all < t < T v , we have a constant d = 0((1 + 7)/e) such that 

f(yn(t))=J2 d " f ~ (X f t)) (yn(t)-x(n(t)r + 0(T-''), 

the constant involved in O only depends on /. We use corollary (|6.ip to get that 

d d 

f{yn(t)) = £ ^/>n(t))(- ^ + 0(T-?)) + O(T^). 

|^|=o 2=0 

Hence we have constants Cx,y,p 1 such that 

f(yn(t)) = £ E w/^%^ + 0(^)- 

|/9|<d,0ez3 z=o 
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Here the O constant only depends on /. Now integrating, we get the result. □ 



Proof of lemma \573[ By using exactly the same proof as in the proof of lemma 14731 
we get that we have constants e x ,yS,i and d! = 0((1 + -f)/e) such that 

f(yu(u)) = J2 E Z^uWfiMu)) + 0{{d' fT-r). 

Integrating on both sides, we get the result. □ 

6.1. Computing c XtVt p t i and ex iVi pj. Recall that d = 0((1 + 7)/e) and we need 
to compute Cxy p,i for all \/3\,l < d. Recall that using (|6.6[) . we have 

\ - fe5 1 ) ,,,„(o)feS, J , I/ (o)4%,x, i/ (o) 

Cx,y,0,l — 2.^1 Qt\ 

P'ez%,\0'\=i 

Here 

r(p - fr) + (tr + s)((p - s)t - f 2 r + g) g 
ni,p u x,yXt) = { TT77 — i — ^2 ) > 



and 



l2 ,fe,^(i) = (-^log(r 2 + (tr + S ) 2 ))^ 



h 3 ,f, 3 , x , y (t) = (tan-^-— 



It is easy to see that we can compute Cx y p i for any Z < d in O(l) time. The 
constant here is a polynomial in d. 

A similar method will go through for computing e x ,y,$,l- 

7. Numerical integration for analytic functions 

Let us state and prove the following simple algorithm for numerical integration 
of a real analytic function, used repeatedly in the paper. In practice this integration 
can be expediated in a number of ways by using a doubly exponential integration 
or Gaussian quadrature technique for numerical integration. 

Proposition 7.1. Let T be any positive number and let f be a real analytic function 
on an open set containing [0, T] . Let I be a fixed integer and R be a positive constant 
such that for any k <G N of \d k f(x)\ < k\R k (l + x l ) for all x in [0,T]. If at any 
given point in [0, T] and n G N and any 7 > 0, we can compute n th derivative of f 
in polynomial (in n ) time, then for any given £,7 > 0, we can compute f(t)dt 
up to an error at most 0(T~ 7 ), using at most 0(T 1+C ) arithmetic operations. The 
constant involved is polynomial in _R, (1 + j)/e and I. 

Proof. The idea is to use a fine grid of T 1+e R cquispaced points and use power 
series expansion around the nearest left grid point to compute /. In particular, let 
us split the integral into integrals over intervals each of size T~ e /R. 
Let 

M = T~ e /R, M 2 = \T/M\ - 1. 

Hence 

pT M 2 pM pT 

(7.1) / f(t)dt = J2 f(xM + t)dt+ f(t)dt. 

JO x=0 J{) JM(M 2 + l) 
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Let us use power series expansion around each xM to compute the value of / at 
a nearby point. Hence we get 

f{xM + t)=Y^d l {f){xM)-. 

1=0 

For \t\ < M and any fixed non negative integer N, we get 



I , „_„ —1-eN 



(7.2) \f{xM + t)-Y,d\f){xM)-\^T l £ T"" 1 « T 

n=0 ' n=N+l 



Substituting this in eauation (17.1|) . we get 

A/2 N M .fc 



/ /(*)*= EE / ^(/)(xM)-dt 

JV /-T-MM2-M fk 

+ Y J d k f{M{M 2 + l)) u dt + E 

k=0 ^° 
M 2 AT „Af .fc 

= EE 9fe (/)(^)/ 

JV pT-MAU-M .fc 

(7.3) +^a fe /(M(M 2 + l)) / t^+£. 

Here, 

|£| « (M 2 + 1)T'- We . 

Notice that any d k f(x) can be computed at any x in polynomial in A; time. Notice 
that the total number of operations needed to compute the sum on the right hand 
side of is 0(M 2 N). We choose N = (1 + i/l) e to § ct thc result - n 

Appendix A. Special functions and lemma 13.21 



As mentioned earlier, the proof of lemma 13.21 given here can also be found, in 
a slight different form in 0. To keep the paper somewhat self contained, we will 
prove it again in this section. 

Let us recall the definitions and some properties of some special functions and 
prove proposition 13.21 



Definition A.l. Given, a, b, c any complex numbers and \z\ < 1, the hypergeomet- 
ric function F(a, b, c, z) is defined by 



1 (all (0) I (n + c)n\ 

There is an analytic continuation of this function to the whole complex plane 
except a branch cut from 1 to infinity. 

We will use the following well known transformation property of Hypergeometric 
functions (see [HI (9.132)]). 
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F(a, b, c, z) = S^r— 4(1 - z)~ a F(a, c-b,l + a-b, 1 



r(6)r(c-a) v ' v ' ' ' 1-z' 

(A.2) + j^ffi ~ S (1 ~ ^"^(6, c ~ o, 1 + & - a, ■ : 



r(a)r(c-&) v ' v ' ' ' i - z y ' 

We also need the following asymptotic expansion for large arguments of the 
Gamma function ([T5J 8.327]) given as 



(A.3) T(z) = z z -ie- z V2ir(l + 0(l/\z\)). 

Hence for z = a + it, a fixed and large |f | we get, 

(A.4) T(a + it) = e l * \t\° )* 1 + O (r 1 )). 

e |t| 

Here, i is positive and the constant in O depends on <r. We use the following result 
about the Besscl functions ([TBI (6.699)]): 

Fact A.l. Let r be any fixed complex number such that \Im(r)\ < \ and T,T' > 0, 
we have 

(A.5) r co,{ Tl t)K ir {t)f T -^dt = 2 ^-V2 T{ i L±E±l n I^±E±l ) 
Jo ll 

ir + iT + | —ir + iT+\ 1 „ 

F( 2-^'2- Tl) - 

Let us use (|A.2[) to get 
(A.6) 

ir + iT+\ -ir + iT+l 1 „ 
F( 2 ' 2 ' 2 ' 1 } 

= r(^)r(-zr)(l + 7f) ~'T~* jr-HT+j jr-jT + j . 1 

T^TjirKl + T*)^ 11 iT-ir+\ -ir - iT + \ _. 1 
T{ lr+l l +h )T{ lr - l l +1 * ) 2 ' 2 ' 1 + 

Now 

F(a,b,c,z) = 1 + 0(|a6z/c|) 

uniformly for 

I | | (ffl + *)(& + /) | . 1 

Lz max — < - 

1 1 ;>o \ c + l)(l + iy ~ 2 

(see [6] and [16]). Hence we can get a positive constant B depending only on 

r such that for any B' > B, both \F( lr+l ^ + 2 , — 2 , ir + 1, 1+ g,-i T i ) — 1| and 

\F( lT !, r+2 , lr 2 ' ~~ + 1' i+b'^t 2 ) ~ 1| arc ^ css than or equal to 0.1. Applying 
this and (|A.6[) to (| A.5[) we get a constant C > B such that for T\ > CT we have 

/>oo 

(A.7) / cos{T 1 t)K ir {t)t T -^dt = D{T ll T){l + r 1 +E{T u T)(l + r 2 )). 
Jo 
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Here |n|,M < 0.1, 

rA R , nrT rp\ 2^r(^±f±i)r(i)rHr)(i + rf) : 

(A.8) D(T X , T) - ir _ iT+ i 

T( 2 -) 



and 



ir+iT+i NtV —ii — iT+± Wl , rr2\ir 



(A.9) £(Tl , T) = .W(-^)r(^-)d + i?) 



r(-ir)r( jr+i 2 T+ ^ )r( ir - i J + g ) 



rfir)r( " iT,+iT+ ^ )r r iT " ir+ ^ ) 
Notice that for a fixed T and real r, the argument of — — — - — . j i . _. ; L i is 

r(- t r)r( ' r+ ' 2 + ^ )r( ' T - ' 2 + ? ) 
fixed. On the other hand (1 + T 1 2 ) lr is a rapidly oscillating function of T\. Hence 
we can choose a suitable C > B such that |1 + n + E(C\T, T)(l + r 2 )| > |. 
Similarly if Re(«r) is non zero, then using asympotics of the T function to (|A.9|) . 
we can see that the behaviour of the absolute value of E(T\,T) is dominated by 
(1 + T 1 2 ) rr , for T\ > T. In other words, we can choose a suitable C > B such that 
|1 + r-i + E(C'T, T)(l + r 2 )\ > \- Moreover, it can be computed in 0(1) time. See 
8.1]. 

Using the asymptotics of T function, it is easy to get a simple bound \D\ ^> T _1 . 
The constant only depends on r. We have thus proved lemma [ 
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